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Probability densities of a forced probe particle in glass: results from mode 
coupling theory and simulations of active microrheology 
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Abstract 

We investigate the displacements of a probe particle inside a glass, when a strong external force is applied 
to the probe (active nonlinear microrhcology) . Calculations within mode coupling theory are presented 
for glasses of hard spheres and compared to Langcvin and Brownian dynamics simulations. Under not 
too strong forces where the probe remains trapped, the probe density distribution becomes anisotropic. 
It is shifted towards the direction of the force, develops an enhanced tail in that direction (signalled 
by a positive skcwncss), and exhibits different variances along and perpendicular to the force direction. 
A simple model of an harmonically trapped probe rationalizes the low force limit, with strong strain 
softening setting in at forces of the order of a few thermal energies per particle radius. 
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1 Introduction 



Many colloidal suspensions are composed of constituents that have typical 
dimensions in the [im range. One is therefore interested in understanding 
the dynamical processes occurring in these systems on mesoscopic length 
scales. A recent experimental technique to obtain such information is 
called active microrheology [l|, Q : one inserts a probe particle, similar in 
size to the host particles, into the system and subjects it to controlled 
external driving, such as a (constant) external force. This is conveniently 
realized in colloids by laser tweezers, magnetic forces, or carefully tailored 
surface-chemistry reactions 0]. Since scaling arguments imply that ther- 
mal fluctuations correspond to forces in the pN range for these colloidal 
systems, one can easily access the strongly nonlinear-response regime in 
active microrheology. If analyzed properly, the nonlinear response of the 
probe reveals a wealth of information on the dynamics of the host fluid; 
in particular if the latter shows highly collective relaxation behavior, such 
as at high densities close to the glass transition Q or close to jamming 
1 • 

In the fluid, a pulled probe eventually reaches a steady-state velocity, 
related to the driving force by a friction coefficient £. Linear response 
demands the velocity to be proportional to the force, hence the friction 
coefficient to be a constant. In active-microrheology experiments on dense 
colloidal hard-sphere model suspensions [fjj], a strong nonlinear variation 
of £ with the external force F ex was found: strong forces drastically reduce 
the friction experienced by the probe. The threshold for this "force thin- 
ning" effect to set in is at forces much larger than the thermal ones, i.e., 
F ex 3> k^T /a, where a is a typical particle radius and k^T is the thermal 
energy. This has been related to the cage effect that dominates the relax- 
ation dynamics at high density: the probe particle is trapped by transient 
cages of nearest neighbors, and the force thinning threshold marks the 
point where force-induced motion overrules the structural relaxation set 
by the decaying cages. 

The limiting case of this scenario is the response of a pulled probe 
inside either a jammed state Q or an ideal glass, i.e., an amorphous solid 
where the embedding particles no longer relax to thermal equilibrium. 
The mode-coupling theory of the glass transition (MCT) 0] explains a ba- 
sic mechanism behind nearest- neighbor cages becoming permanent. The 
theory has recently been extended to force-driven active microrheology in 
the nonlinear-response regime 0]. A remarkable consequence was found: 
in the ideal glass, a force threshold -F e c x exists, which marks a transition, 
from a strictly localized probe in a deformed amorphous solid at low forces, 
to a local melting leading to derealization and long-range motion of the 
probe at large forces. Experiments and computer simulation results in 
the liquid have successfully been interpreted in the theoretical framework 
@, Gil ■ In recent work [II], [lij , emphasis was placed on understanding the 
delocalized state, in terms of friction coefficients and fluctuations around 
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the mean trajectory in terms of mean-squared displacements. 

In this contribution, we focus instead on the localized region, i.e., on a 
probe embedded in a frozen random environment, subject to a force that 
is strong enough to leave the linear response regime, but small enough to 
keep the probe particle localized. This supplements earlier work studying 
single particle motion in crystalline environments [l3j | to the disordered 
case. We present results of MCT describing the probability densities of 
probe positions, and discuss some of their basic statistical properties. We 
compare with results from molecular-dynamics and Brownian-dynamics 
computer simulations. 



2 Methods 
2.1 Theory 

Let us recall the relevant equations from the microscopic MCT for active 
nonlinear microrheology Q . The theory considers the Fourier-transformed 
probe-particle form factor fi at wave vector q. It is finite for a probe par- 
ticle that is trapped in the glass, where it describes the localized probe- 
density distribution at (infinitely) long times. MCT predicts a transition 
from a trapped to a dclocalizcd probe upon increasing the external force, 
or decreasing the particle interactions. The latter can result from decreas- 
ing the host-fluid density, from increasing the temperature of the system, 
or from making the probe particle smaller. Schematically simplified ver- 
sions of MCT have been used to analyze such force-induced derealization 
transitions in detail [l4| , and to describe available data such as friction co- 
efficients, density correlation functions, or mean-squared displacements of 
the probe particle 0, [l(| 12 1. These analyses point out that the idealized 



concept of a forced but still localized probe in a non-relaxing amorphous 
host solid is very useful in understanding the qualitative nonlinearities 
observed for probes pulled through dense, complex fluids. 
The real-space representation of the form factor fi, 



describes the stationary probability for the center of the probe particle 
to be in an infinitesimal volume element around r, provided that it was 
at the origin initially (when the constant external force F ex was switched 
on). Because of particle conservation, f s (r) is normalized to unity. In a 
quiescent amorphous solid, spatial isotropy is assumed, which is reflected 
in the fact that /| depends on the wave vector only through its modulus, 

\o\- 

Spatial symmetry is broken by the external force: the spherical sym- 
metry is reduced to a rotational one around the force axis (taken to be 
the z axis in the following), and the distribution in the direction of the 
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force becomes qualitatively different. In particular, it allows for a non- 
zero mean and, more generally, non-zero odd moments of the distribution 
- something forbidden by symmetry in the quiescent case. For later ref- 
erence, recall the standard definitions of the first few moments of f s (f): 
the lowest nontrivial moments are the mean (z) along the force direction, 
the variances (Az 2 ) and (Ax 2 ) in direction of and perpendicular to the 
force, and the skewness of the distribution along the force axis, 71: 

(2a) 
(2b) 
(2c) 

dr (z-{z)f r(f) (2d) 



(z) = 


J dfzf s (r) 


(Az 2 ) = 


J dr{z-(z)f f s (r) 


(Ax 2 ) = 


J drx 2 f(r) 




(Az 3 ) 1 



(Az 2 )3/2 (Az 2 ) 3 / 2 

MCT decomposes the internal forces acting on the probe into density 
fluctuations of the probe and its surrounding host liquid, which we take 
to consist of a colloidal dispersion dissolved in a (featureless) solvent. 
Because the particle interactions are translational invariant, MCT works 
with plane- wave decompositions of these fluctuations, and only at the end 
Fourier back-transforms to real space are taken as in Eq. ([T]). The form 
factor fi appears as the long-time limit of probe-density autocorrelation 
functions, that do not decay to zero for a localized particle (corresponding 
to the emergence of an elastic contribution in the scattering spectrum). 
Quite generically, one derives 

m s -.\f, f s ] 

The cage effect appears in the retarded friction and its history-dependence. 
It is captured in a generalized friction kernel ml for which exact formal 
expressions are known that, however, can only be evaluated with ap- 
proximations. MCT approximates the host-probe force fluctuations using 
effective interactions provided by the Ornstein-Zernike direct correlation 
functions. As a consequence, only wavevector-dependent density fluctu- 
ations need to be calculated within the theory. The direct correlation 
functions here are assumed to be known from some separate theory, e.g., 
density functional or liquid state theory (IBj . In the idealized glass, the 
structure of the host system is then fully encapsulated in f q , the isotropic 
collective nonergodicity parameter (Debye- Waller factor, in the case of 
a monodisperse one-component host). The f q are real- valued and well 
known from quiescent MCT calculations Q . 

The wavevector-dependent probe form factor /|, the central object of 
our calculations, then is determined by the MCT closure to Eq. ([3]), given 
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through the mode-coupling functional 



Here, p = q — k, and p = \p\ have been introduced for abbreviation, 
and the host-liquid particle density is denoted by n. The term cot - = 

[q — iF cx /kBT} ■ p embodies instantaneous forces. S s (p) and S(p) are the 
equilibrium static structure factors given by, respectively, the probe-host 
and the host-liquid interactions, and connected to the mentioned direct 
correlation functions. The latter structure factor also is the only input 
into the calculation of the collective noncrgodicity parameter f q within 
MCT, so that all input to determine fi is now specified. 

While for F cx = the well-studied quiescent MCT follows, finite ex- 
ternal forces F cx affect by Newton's laws the momentum balance and 
thus cause a shift in the wavevector q in w^p-. Because of the displace- 
ment of the probe from the origin, its Fourier-transformed distribution 
/| is complex. Realness of the real space distribution, however, follows 
from fi* = /£.-, which is a consequence from Eqs.([3]) to (|4]), at least 
for small and intermediate forces. There arc indications from numerical 
solutions of the microscopic MCT equations that at larger forces, a bi- 
furcation takes place to a solution that violates this requirement |16j |. In 
the previous analyses involving schematic models, where no detailed spa- 
tial information is kept, this singularity was regularized. Here, we only 
consider external forces that are smaller than the bifurcation threshold. 

Choosing the interactions between particles determines the equilibrium 
structure information required in Eq. (j4|), and thus completely specifics 
the solutions of the MCT equations for active microrheology. We will 
discuss a system of monodispcrsc hard spheres with radius a. The di- 
mensionless density is given in terms of the packing fraction, ip = 4^na 3 , 
which quantifies the volume taken by the spheres relative to the total 
volume. It is the only control parameter determining the structure of 
the host system. For simplicity, we will also consider a probe particle 
that is identical in size to the surrounding spheres. The Pcrcus-Yevick 
approximation [lBj of the direct correlation functions for hard spheres is 
used to approximate the MCT coupling coefficients: this liquid-state the- 
ory is well established, analytically solvable, and remains reasonable up 
to high densities. Within this model, MCT predicts a glass transition to 
occur at a critical packing fraction ip c « 0.516. We consider the window 
fc < ip < 0.54 in the following. 

Equations and are solved on a discrete lattice of 128 x 129 
points in the magnitude-angle plane, (|g|, ■&), where ■& is the angle formed 
between q and F cx . Angular integration runs over the interval [0,7r] for 
the angle fl. As fi is rotationally invariant around the axis of the force it 
does not depend on the azimuth </>. Still integration of this angle over the 
interval [0, 2tt] has to be performed when evaluating the mode-coupling 
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functional. Wave-vector modulus integration starts at 2ag m j n = 0.001 up 
to a cutoff 2aq max = 45, 50, 60, or 65 (for ip = 0.516, 0.52, 0.53, and 
0.54); at larger q, both f q and /| have sufficiently decayed to zero. To 
allow for sufficiently smooth interpolation of f q , which is pre-determined 
by equations analogous to Eqs. ^ and ((H), the latter was calculated on a 
grid of 2048 points in \q\ (angular integration being performed analytically 
in this case). We found no qualitative dependence on the discretization, 
checked by some calculations for finer and coarser grids. However, the 
inverse Fourier transform required by Eq. (fTJ) is very sensitive to cutoff 
effects, and for the grid used here causes tiny rapid oscillations in f s (f) 
that cause artefacts when calculating moments of the distribution. This 
is remedied by setting |/ s (r)| to zero below a threshold of 2.5 10~ 4 . 



2.2 Simulations 

To conduct first qualitative tests of the predicted features of the density 
distributions, we also performed computer simulations of active force- 
driven microrheology. Two model systems were used: strongly damped 
molecular dynamics for a three-dimensional system of slightly soft spheres 
based on a Langevin equation, and two-dimensional Brownian dynamics 
simulations for hard spheres. 

The Langevin dynamics simulation has been used extensively before, 



both in discussing the equilibrium approach to the glass transition jl7l.ll 



squill brn 

and for microrheology [3, [l(| ; implementation details are given in these 
publications. The 3d system consists of 1000 purely repulsive particles 
with a power-law interaction, V(r) ~ r -36 , and interaction radii that 
are randomly drawn from a flat distribution of half-width 10% to avoid 
crystallization. Periodic boundary conditions are used, and the box is 
elongated by a factor of 8 in the force direction to reduce finite-size effects. 
The glass-transition packing fraction of this system is known to be tp c ~ 
0.595. Here we show simulations sweeping the external force at a fixed 
density of cp = 0.57, accessing the (still transiently) frozen-in component of 
the probe density distribution at an intermediate time t w 25y // ma 2 jh&T 
(m is the particle mass). For this time, the transient self- intermediate 
scattering functions at qa « 3 are close to their plateau value 0, [Toj j . 
Averaging over 300 simulations each consisting of about 40 independent 
driving cycles was used to improve statistics. Simulations have also been 
performed at lower densities, and at tp = 0.62; for the latter, the system 
can no longer be equilibrated in the available computing time, but has 
been aged for t w = 25000 ^JTtmFJk^f . Still, with this choice, considerable 
aging effects cannot be ruled out for F ox < 35 k^T/a, making it preferable 
to obtain the density distributions from the liquid-state simulation as 
mentioned above. 

The Brownian dynamics simulation considers an equimolar binary 
mixture of hard discs in d = 2 dimensions, with radius ratio 1.4. A hybrid 
Monte-Carlo scheme based on event-driven dynamics fljj is used, where 
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Figure 1: Contour plot of the probability density of a forced probe-sphere 
in a hard sphere glass at the glass transition packing fraction, ip c = 0.516, 
as calculated in MCT. The probe-particle and host particles are equal 
sized (radius a), and the force is F ex = lQe z k-Q,T/a. A cut f s (x,z,y = 0) 
is shown. 



a thermostat is introduced which at every integer-multiple of a 'Brow- 
nian time' tb redraws all particle velocities from an identical Gaussian 
distribution. This mimics Brownian motion after many tb- Inspired by 
Carpen and Brady [2(J, the external force on the probe-particle is imple- 
mented by a shift of its Gaussian velocity distribution by the value F ex / Co 
in the fixed force direction (Co is the short time friction coefficient). The 
glass-transition dynamics of this system has been analyzed in detail [2l| . 
establishing tp c ps 0.795. For the reasons explained above, we chose a den- 
sity <p = 0.78 and measured the probe-particle distribution function at 
intermediate times. As a probe particle, a bigger disc is chosen randomly. 



3 Results 

Figure Q] shows a two-dimensional cut, f s (x,z,y = 0), through the prob- 
ability density of the forced probe in a hard sphere glass as calculated 
within MCT. The packing fraction is the critical one, <p c = 0.516, where 
quiescent MCT locates the formation of an (idealized) glass. Probe and 
host fluid particles are taken to be identical, which leads to a strong cou- 
pling of the probe-particle to the collective fluid density. For vanishing 
external force, the probe is localized inside nearest-neighbor cages, and 
due to spatial isotropy its probability density is spherically symmetric and 
centered around the origin. The width of the distribution is of the order 
of a tenth of the average particle separation (that is, some fraction of 
the particle diameter 2a), as was first discussed by Lindemann in crystals 



113,123. 

For the large force F ex = 10 UbT /a shown in Fig. [IJ the distribution is 
shifted to the wall of the cage, and a slight anisotropy between the direc- 
tions aligned with and opposite to the external force becomes noticeable. 

To identify the anisotropic effects induced by the external force more 
clearly, it is convenient to discuss the one-dimensional marginal proba- 
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Figure 2: Stationary density distributions of the probe-particle along the 
direction of the external force as calculated by MCT, and for forces as 
labeled in units of k^T / a. Probe and fluid particles have radius a. The 
upper panel is at the packing fraction of the MCT glass transition, ip c = 
0.516, the lower one in a glass at (p = 0.53 (note the different horizontal 
scales). The insets show the corresponding probability densities along 
the ^-direction perpendicular to the force (with again different horizontal 
scales). 
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bility densities obtained from integrating out two spatial directions. Of 
particular interest is f s {z), the probability distribution along the force 
axis, shown in Fig. [5] for forces starting from zero up to the value F cx = 
10 k^T /a shown in Fig. [TJ in steps of 2.5 ksT/a. In the insets of Fig. [2J 
we also show f s (x) = f s (y)-, the marginal distributions along directions 
perpendicular to the force. We show results for both the glass-transition 
density, and a state point deeper in the glass, if = 0.53, where the probe 
remains localized for a larger range of forces. The qualitative effects of the 
external force are the same irrespective of density (as long as one stays 
inside the glass). Quantitatively, the distributions for the higher density 
are narrower, reflecting the fact that particles are more tightly localized 
in a more densely packed system; note the different horizontal scales in 
the two panels of Fig. [2j 

With increasing force, the half-width at half maximum of the distri- 
bution increases, starting from Lindemann's value at zero force. This 
indicates that cages are continuously widened by applying an increasing 
external force, in all spatial directions. The center of the distribution 
along the force axis, f s (z), shifts in the direction of the force, as is in- 
tuitive. Spatial symmetry demands that f s {x) at the same time remains 
centered on zero, as verified from the figure. 

In f(z), a characteristic asymmetry develops, with an increased wing 
extending in the direction of the force stronger than in the direction op- 
posite to it. Interestingly, also the distributions characterizing the probe 
probability densities in directions perpendicular to the force show ap- 
preciable wings at large forces, extending to distances much larger than 
Lindemann's length. 

Let us remark on one other feature that can be identified for the distri- 
butions shown in Fig. [2j For some state points, slightly negative dips can 
be identified in the MCT calculations. This is clearly unphysical. One has 
to keep in mind that MCT is formulated in wave-vector space, where it 
is designed to guarantee a number of exact properties of correlation func- 
tions in equilibrium Since the theory involves approximations, this 
guarantee however does not carry over to the real-space formulation ob- 
tained from Fourier-transforming the MCT results. In fact, small negative 
contributions in the predicted probability densities (van Hove functions) 



of a quiescent hard-sphere glass arc known |18| : they are thus not a con- 
sequence of extending the theory to the nonlinear response. Note that 
even the equilibrium structure approximation we use as input for our cal- 
culations, viz. the Pcrcus-Yevick approximation for hard spheres, fails to 
satisfy exact positivity of the real-spa ce p air distribution functions (albeit 
at higher densities than we consider) [15j. As Fig.[2]shows, the unphysical 
negative dips become more pronounced for higher densities, but become 
weaker and finally (almost) disappear with increasing force. They remain 
small for all densities ip < 0.54. In the following, we will concentrate on 
the upper portions of the distribution functions where we expect MCT to 
make sensible predictions. 
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ex 1 B ' 

Figure 3: Mean value (z) of the long-time limit of the probe-particle 
displacement along the force direction as function of the external force 
F cx . Curves for different glass packing fractions ip are shown as labeled. 
At the higher two ip, corrected (filled squares) and uncorrected (open 
squares) (z) differ slightly; see text for details. 

From the overall shape of the distributions shown in Fig. [2] it is ap- 
parent that their first few moments allow to quantify the mentioned qual- 
itative effects. We now turn to a discussion specifically of the mean, 
variance, and skewness of f B (z) and f s (x). In calculating these quanti- 
ties, the unphysical negative portions have to be treated with care. In 
particular high-order moments tend to emphasize the small errors iden- 
tified in Fig. [2] We thus introduce an ad-hoc correction by setting the 
negative dips to zero when calculating moments (the distribution is re- 
normalized afterwards). Values obtained from both the uncorrected and 
the corrected distributions are shown in the following figures, to allow an 
assessment of our ad-hoc correction. 

Figure |3] shows the mean values (z) of the probe-displacements in re- 
sponse to the external force F ox for a number of glass states. The curves 
are continued up to the values of the external force where the mentioned 
bifurcation in the MCT equations takes place. Beyond this force, which 
depends on density, we consider the microscopic MCT results not trust- 
worthy. With increasing density, the shift of (z) with F cx becomes smaller 
as the elastic strength of the glass also increases. Initially, (z) grows lin- 
early with F cx : this may be called Hookean behavior, and is expected 
within the general approach of linear-response theory. We discuss this in 
more detail below, in connection with Fig. 1101 Approaching the bifurca- 
tion point, a highly nonlinear dependence of (z) on F ex is seen. Since (z) 
grows stronger than linearly, one may call this "strain softening" . (Reas- 
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F ex /(k B T/a) 

Figure 4: Variances (Az 2 ) (squares) and (Ax 2 ) (triangles) of the station- 
ary probe-particle density distribution along, respectively, perpendicular 
to the force direction as function of the external force F ox . Curves for dif- 
ferent glass packing fractions cp are shown as labeled. At the higher two ip 
the corrected (filled squares) and uncorrected (open squares) (Az 2 )-data 
differ slightly; see text for details. 



suringly, the effect of the mentioned ad-hoc correction is small and only 
becomes visible for the highest densities considered here.) 

Figure 2] shows the variances defined in Eq. ^ as functions of the 
strength of the external force for the four glass densities used above. 
Again, corrected and uncorrected values agree quite well except for small 
differences in (Az 2 ) at the highest densities. At small F eK , the probe- 
particle response remains nearly isotropic, as required by linear response. 
Both (Az 2 ) and {Ax 2 } approach for F ex — > values given by the localiza- 
tion length of the quiescent glass, r s . This quantity has been evaluated 
before [8| and has been connected to the Lindemann criterion already men- 
tioned: the localization length r s of single particles inside a crystal cannot 
exceed some fraction of the lattice constant, as l ong as the crystal is sup- 



posed to be stable against thermal fluctuations [22J. MCT transfers this 
concept to a disordered solid, where the average particle separation plays 
the role of the relevant length scale. The value (Az 2 ) = 2r 2 « 0.04a 2 
at the transition density tp c has already been verified by dynamic light 
scattering experiments on hard-sphere- model colloidal supsensions [231.124 1 
and by simulations [13] • In these cases, the time-dependent mean-squared 
displacement (Sr 2 (t)) was studied, whose long-time limit in the glass is 
{5r 2 (t)) = 3(Az 2 ) assuming spatial isotropy. 

For increasing force F CK , Fig. [4] confirms the qualitative judgement of 
Fig. [2j the probe-particle distributions become wider in all directions, 
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F ex /(k B T/a) 

Figure 5: Skcwncss 71 of the probe-particle density distribution along the 
force direction as function of the external force F cx . Curves for different 
glass packing fractions ip are shown as labeled. At all (p except for the 
lowest one, corrected (filled squares) and uncorrected (open squares) 71 
differ, because the latter arc affected by the (small) negative tails in f s (z) 
predicted by MCT. Thus, we judge the corrected ones more reliable; see 
text for details. 

revealing a rather steep increase as the bifurcation force is approached 
from below. Outside the linear-response regime, the response also be- 
comes anisotropic. The distribution perpendicular to the force direction 
becomes wider more quickly than along the force. We attribute this effect 
to the wide tails seen in the insets of Fig. [2] 

We next discuss the skewness of the distributions, Fig. [5] Here, the 
small negative dips in the MCT approximation to the real-space distri- 
bution functions are amplified so that the uncorrected values of 71 be- 
come strongly negative for large densities and strong forces. The ad-hoc 
corrected values remain small at small forces and then are positive and 
increase appreciably for larger F ex ; since we expect MCT to work best for 
the central features of the distribution we consider these latter values to 
be more accurate. Generically, the marginal distribution f s (z) is hence 
right-skewed, i.e., its mass is concentrated to values lower than the mean, 
indicating that a few realizations correspond to large excursions in the 
direction of the force. This is most prominent close to the bifurcation 
point, where again the most steep increase is seen. 
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Figure 6: Linear displacement of the forced tracer particle in direction of 
the force (given in k^TJa) for the density (p = 0.57 from Langevin dy- 
namics simulations. The intermediate time t = 25y 'ma 2 /k-gT at which 
probe-particle distributions were sampled is indicated by the dashed ver- 
tical line. 



4 Discussion 

4.1 Preliminary comparison with simulations 

In order to compare with computer simulation, we face the problem that 
the idealized glass state considered by MCT cannot be obtained in the 
simulation: at high densities, the system can no longer be equilibrated in 
the available computer time, and aging effects prevail that are outside the 
scope of our discussion. We hence choose a density close to the MCT glass 
transition, where the random host structure still relaxes in the simulation 
time window. Before this final relaxation sets in, the intermediate nearly- 
frozen structure can be probed. 

To demonstrate this, Fig. [6] shows the linear displacement (z)(t) of 
the probe particle as a function of time since switching on the external 
force, obtained from the Langcvin-dynamics simulation. One recognizes 
that these functions grow linearly with time at large t, indicating fluid-like 
steady-state motion corresponding to a finite friction coefficient £. This 
steady state is reached on a time scale governed by the host-fluid relax- 
ation at small forces. Before that, (z)(t) exhibits a shoulder as a function 
of log t, indicative of transient localization. We pick an intermediate time 
t = 25^ma 2 /ksT (indicated by a dashed line in Fig. [6]) where probe- 
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Figure 7: Probability density of a forced probe in a polydisperse hard 
sphere fluid at ip = 0.57 from Langevin dynamics simulations. From top 
to bottom: F ex /(kBT/a) = 10, 20, and 30; here the force is along the 
x-axis. The colour code for the unnormalized distribution is given in the 
upper panel. The time since switching-on the force is chosen such that 
the localized density distribution is close to its frozen-in form at the glass 
transition. 

particle distributions arc sampled. At this time, the probe particle has 
fully explored its cage, but is still sufficiently close to its transient localized 
state for a range of low forces, roughly F cx < 2Qk-oT/a. For the higher 
forces, the time when the probe is pulled out-off the cage moves to shorter 
values, and the probability distribution sampled at t — 2h\J 'ma? jk^T be- 
comes characteristic for the delocalized state. For F ox > 2QkBT/a at the 
given packing fraction (ip = 0.57) the simulation results for f 3 thus will 
contain different phenomena than discussed from the theory. 

Figure [7] shows the probe-particle distribution f s (r) sampled at that 
time, for various external forces. As in Fig. [5J a cut in a plane containing 
F ex has been chosen, and rotational symmetry around that axis has been 
used to improve statistics. The panels correspond to increasing external 
force. One notices the nearly intact spherical symmetry for the lowest 
force, F cx = 10 k&T/a, which is broken for the larger forces. Qualitatively, 
one recognizes the features discussed above from the MCT predictions: 
an increasing shift of the center of the distribution in the direction of the 
applied force, and an increase in the variance in all spatial directions. 
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Figure 8: Langevin simulation results for the skewness 71 of the proba- 
bility distribution along the force-direction as a function of the force, for 
different densities as labeled. Forces are given in units of k^Tja. 



The heterogeneity of the probe dynamics is parametrized by the skew- 
ness of these distributions in Fig. as defined in Eq. (|2dj) . Increasing the 
external force from zero, the skewness increases, before it reaches a max- 
imum at intermediate forces and finally decays to zero again. From the 
above discussion of Fig. [S] and a comparison with the mean probe velocity 
H, we associate this decay with the eventual derealization of the probe; 
this effect is not contained for the small forces discussed within MCT 
above. There is qualitative agreement between theory and simulation for 
the initial increase: The skewness is positive, indicating that the probe 
undergoes rare excursions far in the direction of the force, i.e., the density 
distribution develops a tail there. While in the fluid, the maximum of 
71 lies at the same force, in the glass the increase shifts to larger forces. 
This is compatible with the MCT scenario of a critical bifurcation force 
whose value in the fluid close to the glass transition is given by the one 
at the transition. Values of 71 « 1 at forces around F ex « lQk-BT/a at 
the glass transition density are found in MCT; this is in semi-quantitative 
agreement with the simulations. 

Also the two-dimensional Brownian dynamics simulations of active 
microrheology, broadly give results in agreement with the theoretical pre- 
dictions. Figure [9] shows a contour plot of the density distribution of a 
probe (one of the larger species) trapped in an equimolar mixture of disks 
with radii-ratio a,b/a s = 1.4. The density is close to the glass transition 
one, and a large force F CK = 42/cbT /a is applied. Here a denotes the aver- 
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Figure 9: Probability density of a forced probe in a binary hard disk 
mixture (radius ratio ab/a s = 1.4) in d = 2-dimensions close to the glass 
transition density from Brownian dynamics simulations. The probe has 
been chosen among the bigger disks, the applied force is F cx = ^Ik^T ja 
in z-direction and a = 0.6 denotes the average particle radius. The colour 
code gives the probability. The time since switching-on the force is chosen 
such, that the localized density distribution is close to its frozcn-in form 
at the glass transition. 

age particle radius a = (a& + a s )/2 . Estimates of the low-order moments, 
(z) w 0.14a, (Az 2 ) w 0.32a 2 , (Ax 2 ) « 0.24a 2 , and 71 w 2.51 lead to 
values which can be taken as extrapolations of the theoretical trends to 
such strong forces; the small average displacement is noteworthy though, 
and presumably due to the lower spatial dimension. The asymmetry of 
the density distribution along the force direction is apparent already from 
the contour plot. 

4.2 Small force model of an harmonically trapped 
probe 

To place the findings reported above in theoretical context, we briefly dis- 
cuss the linear-response theory of a forced particle in a glass. A prototyp- 
ical model system is that of a Brownian particle in a harmonic trap. Con- 
sidering for simplicity only the force direction, there holds the Langcvin 
equation 

C«(t) = -k (u(t) - uo) + F cx + f(t) . (5) 

Here, u(t) is the displacement variable of the particle, and k the spring 
constant of the harmonic potential. Random forces / are assumed to 
be Gaussian white noise obeying the fluctuation-dissipation theorem con- 
nected to the friction coefficient £. The minimum position of the potential 
is given by uq. The distance from this minimum, Su = u(t) — uq immedi- 
ately gives the linear shift and the equal-time variance of the displacement 
at infinite times, 
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Figure 10: Effective force constant k (circles) determined from the probe- 
particle displacement under force for different glass packing fractions vp 
as listed in the legend. Testing the linear response relation Eq. ([6|), the 
displacement variances along k^,T / (8u 2 ) gl = 2k^,T / '(Az 2 ) (triangles) and 
fceT '/ {Su 2 ,)^ = 2/cbT / (Ax 2 ) (squares) perpendicular to the force are com- 
pared to k = F cx /(8u) = F ex /(z). For small forces, where (Su 2 ) gl — > r 2 s 
isotropically, linear response and the equipartition theorem hold. Shown 
are only the corrected values not affected by the (small) negative tails in 
f s (z) predicted by MCT; see text for details. 
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Yet, the glass state is not an equilibrium state, but a metastable state, 
where ergodicity has been broken. The nonergodicity parameters f q of 
MCT quantify this, and a finite localization length r s = J (8z 2 (t — >• oo))/2 
is the simplest example of a nonergodic quantity in the probe dynamics. 
Importantly, the localization length is measured in the limit of infinite 
times, where the maximal excursion Sz(t) = z{t) — z(0) of the tracer 
can be obtained. The question arises how to define ergodic, close to 
(metastable) equilibrium fluctuation functions in the glass state, which 
are assumed in the simple harmonic oscillator model. Recently, in the 
context of measuring strain fluctuations in colloidal glass by microscopy, 
this problem has been solved [25[. The stationary (quasi-equilibrium) 
displacement auto-correlation function in the glass state is defined from 
the mean squared displacement defined in the fluid state by 

(5u(t)5u(0))* = \ ((Az 2 ) - <(*(*) - z(0)) 2 )) . (7) 

Here averages on the right are performed using the Gibbs-Boltzmann dis- 
tribution of the fluid state, while the average on the left hand side is 
within the compartment in phase space corresponding to the nonergodic 
glass state. The function (5u(t)£u(0)) gl defined on the left hand side obeys 
all the properties required of an auto-correlation function, at least if the 
approximate MCT equations of motion are taken as basis [25| . Clearly and 
as expected, Eq. (0 gives a linear relation between the displacement u(t) 
of the probe in the glass and the position z(t), whose long time distribution 
functions were the topic of our MCT calculations: u(t) = z(t). Somewhat 
more surprisingly, the second moment of the displacement taken at equal 
time follows to be half as big as the variance of the density distribution 
remaining at long times: 

(SuY = ((6u(t = 0))Y = \ <Az 2 } = \ ((z(t -+«>)- z(0)) 2 ) . (8) 

Figure [TU] tests the thus obtained simple model of the motion of the 
trapped probe in glass. The values of the effective spring constant, k = 
F cx /{z), are compared to the appropriately inverted mean squared dis- 
placements (Su 2 ) sl , namely k^T / {Su 2 ) gl , and quantitative agreement is 
found. This supports the above definition of a displacement fluctuation 
function in the nonergodic glass state. Its behavior at small forces reduces 
to the simple model of an equilibrium harmonic oscillator. At larger ex- 
ternal forces, force thinning sets in, and the first and second displacement 
moments behave qualitatively similar, yet quantitatively differ. (In order 
to eliminate spurious non-monotonic behavior arising from the negative 
probabilities MCT would predict, only the corrected moments, where the 
negative probabilities were set to zero, are shown in Fig. 1101 ) 

The linear-response spring constant evaluated by MCT close to the 
glass transition is k w bOk^T/a 2 , which agrees quite well with a re- 
cent estimate k « 80 k^T/a 2 obtained employing a parabolic trap in the 
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Langevin simulations [261 ] . Following basic elasticity theory, it can be 



related to a (local) Young's modulus E of the glass. Another relevant 
length enters this conversion (the cross-section over which the force is ap- 
plied divided by the original length of the strained object). Assuming this 
to be roughly a again, we get E m i cro ps bOk^T/a?. In the linear macro- 
scopic rheology of glasses, MCT predicts both the longitudinal and the 
shear modulus, M and G; for a hard-sphere system at the glass transition, 
G ps 20fc B T/(2a) 3 and M = M + 5AI ps 130 fc B T/(2a) 3 , where SM ps 
55kBT/(2a) 3 is the nonergodic contribution [27j . From the standard re- 
lation between the elastic moduli, E = (G/M)(l - (4/3)(G/M))/(l - 
(G/M)), one gets E ps 56 fceT/a 3 , which compares reasonably well with 
the microscopic spring constant. 



5 Conclusions 

We investigated the nonlinear response of a single particle subject to a 
strong external force in an amorphous solid. Calculations based on the 
mode-coupling theory for active nonlinear microrheology allow to access 
the density distribution functions describing the average probe position in 
the nonequilibrium system. These distribution functions become asym- 
metric and loose the spherical symmetry present in the spatially homoge- 
neous equilibrium state: quite intuitively, the center of the distribution is 
shifted in the direction of the applied force. In addition, two nontrivial ef- 
fects are predicted: the distribution is widened in all spatial directions, i.e., 
its localization also perpendicular to the applied force becomes less strong 
when the force is increased. Further, one identifies strong tails in the dis- 
tribution function, and in particular a pronounced tail extending along 
the force from the initial positions, as manifested in a positive skewness 
parameter. This tail becomes stronger as the force approaches a critical 
bifurcation threshold. It quantifies dynamical heterogeneities that arise 
because occasionally particle trajectories follow realizations with large 
excursions. In fact, it is seen in the simulation that the probe-particle 
dynamics close to derealization becomes intermittent. The connection of 
the strong tails and the bifurcation force found in the MCT equations is 
currently under investigation. 

For small forces, linear response theory allows to understand the re- 
sults if one takes into account the broken ergodicity of the glass. This 
proper definition of quasi-equilibrium displacement correlation functions 
in glass accounts for an otherwise puzzling factor 2 connecting the variance 
of the probe-distribution functions to the harmonic spring constant char- 
acterizing the local rigidity of the nearest-neighbor cage. Semi-quantitatively, 
the local measure of rigidity can be connected to the macroscopic elastic 
constants of the glass. 
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